clear
clc
inputName = '1.jpg';
X=imread(inputName);
Prtemp=X(:,:,1);
Pgtemp=X(:,:,2);
Pbtemp=X(:,:,3);

Pb=double(Pbtemp(:));
Pr=double(Prtemp(:));
Pg=double(Pgtemp(:));
P=[Pr,Pg,Pb];

[Prow Pcol num] = size(P);
modify_input = zeros(Prow,Pcol);
W = ones(1,Prow);
% W = W + mean(P')/255;
for(i=1:Prow)
    modify_input(i,:) = double(P(i,:))*W(i);
end

[V sigma] = eig(modify_input'*P);
S = sqrt(max(sigma));
SV = (S*V(:,3))';
%B=V(:,3);
B = 255*(1/max(V(:,3)))*(V(:,3));
Br = B(1);
Bg = B(2);
Bb = B(3);
Pr = (1/Br) * Pr;
Pg = (1/Bg) * Pg;
Pb = (1/Bb) * Pb;
[m n] = size(Pr);
T = zeros(m,1);
T = min(min(Pr,Pg),Pb);
[row col num] = size(X);
R =double(T) * B';
cover_energy = sum(sum(R)) / sum(sum(P));

% %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% % scaling the matrix
% %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% if (cover_energy < 0.5)
%     R =0.99/cover_energy*R;
%     red = sum(sum(R(:,1))) / sum(sum(P(:,1)));
%     green = sum(sum(R(:,2))) / sum(sum(P(:,2)));
%     blue = sum(sum(R(:,3))) / sum(sum(P(:,3)));
%     if (red>1)
%         R = (1/red)*R;
%     end
%     if (green >1)
%         R = (1/green)*R;
%     end
%     if (blue>1)
%         R = (1/blue)*R;
%     end
% end
% cover_energy = sum(sum(R)) / sum(sum(P));
% red = sum(sum(R(:,1))) / sum(sum(P));
% green = sum(sum(R(:,2))) / sum(sum(P));
% blue = sum(sum(R(:,3))) / sum(sum(P));
% %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%


Rr = R(:,1)';
Rg = R(:,2)';
Rb = R(:,3)';
RP(:,:,1) = reshape(Rr,row,col);
RP(:,:,2) = reshape(Rg,row,col);
RP(:,:,3) = reshape(Rb,row,col);

imshow(uint8(RP));
ti = '1stField_';
outputName = [ti inputName];
imwrite(uint8(RP), outputName);

Rrc = double(X(:,:,1)) - RP(:,:,1);
Rgc = double(X(:,:,2)) - RP(:,:,2);
Rbc = double(X(:,:,3)) - RP(:,:,3);

Rc2 = zeros(row,col,3);
Rc2(:,:,1) = Rrc;
Rc2N = '2ndField_';
outputName = [Rc2N inputName];
imwrite(uint8(Rc2), outputName);
second_energy = sum(sum(sum(Rc2))) / sum(sum(P));

Rc3 = zeros(row,col,3);
Rc3(:,:,2) = Rgc;
Rc3N = '3rdField_';
outputName = [Rc3N inputName];
imwrite(uint8(Rc3), outputName);
third_energy = sum(sum(sum(Rc3))) / sum(sum(P));

Rc4 = zeros(row,col,3);
Rc4(:,:,3) = Rbc;
Rc4N = '4thField_';
outputName = [Rc4N inputName];
imwrite(uint8(Rc4), outputName);
fourth_energy = sum(sum(sum(Rc4))) / sum(sum(P));

% RP(:,:,1) = RP(:,:,1) + Rrc;
% RP(:,:,2) = RP(:,:,2) + Rgc;
% RP(:,:,3) = RP(:,:,3) + Rbc;
% Rec = 'Reconstruct_';
% outputName = [Rec inputName];
% imwrite(uint8(RP), outputName);
% Reconstruct_energy = sum(sum(sum(RP))) / sum(sum(P))
fprintf('1. cover_energy = %.2f %%\n',cover_energy*100);
% fprintf('   - red_energy = %.2f %%\n',red*100);
% fprintf('   - green_energy = %.2f %%\n',green*100);
% fprintf('   - blue_energy = %.2f %%\n',blue*100);
fprintf('2. second_energy = %.2f %%\n',second_energy*100);
fprintf('3. third_energy = %.2f %%\n',third_energy*100);
fprintf('4. fourth_energy = %.2f %%\n',fourth_energy*100);
%fid = fopen('data.txt', 'A+');
%fprintf('image name    row    col    main field%%    R%%    G%%    B%%    ')
%fprintf(fid,'%s|    %d|    %d|    %.2f%%|    %.2f%%|    %.2f%%|    %.2f%%    \n',inputName,row,col,cover_energy*100,red*100,green*100,blue*100); 
%fclose(fid);
return 